Robustness of One-Dimensional Photonic Bandgaps Under Random Variations of 

Geometrical Parameters 



H. Sami Sozuer and Koray Sevim 
Izmir Institute of Technology, Department of Physics, Gulbahce Koyu, Urla, Izmir, TURKE^ 

The supercell method is used to study the variation of the photonic bandgaps in one-dimensional 
photonic crystals under random perturbations to thicknesses of the layers. The results of both plane 
wave and analytical band structure and density of states calculations are presented along with the 
transmission coefficient as the level of randomness and the supercell size is increased. It is found 
that higher bandgaps disappear first as the randomness is gradually increased. The lowest bandgap 
is found to persist up to a randomness level of 55 percent. 



I. INTRODUCTION 

Since the pioneering work of E. Yablonovitch pj and 
S. John 0, research on photonic crystals(PCs) has en- 
joyed a nearly exponential increase. The manufacture of 
PCs at the optical regime has become a reality Q. Man- 
ufacturing brings with it the practical reality of random 
errors introduced during the manufacturing process and 
it is the effect of these random errors on the desirable fea- 
tures of PCs, namely photonic band gaps, that we wish 
to address in this paper. 

Bandgaps in PCs depend on two crucial properties: an 
infinite and perfect translational symmetry. Clearly, in 
real life no crystal is infinite in size or perfectly periodic. 
When randomness is introduced in the geometry of the 
PC, one quantity of interest is the size of the bandgaps 
as the level of randomness is increased, and whether the 
bandgaps of the bulk perfect PC will survive the random- 
ness. Same considerations apply for a finite PC. In fact, 
even for a perfect but finite PC, one needs to give up 
the notion of a bandgap and has to be content with se- 
vere depressions in transmittance instead. In this paper, 
we will consider both finite imperfect PCs by examining 
the dependence of their transmittance on randomness, 
and bulk imperfect PCs by determining their density of 
states (DoS) under varying degrees of randomness, using 
the supercell method. 

Although much has been done HUE! [13 

regard- 
ing imperfect two- and three-dimensional PCs, we feel 
that a study of the problem for one-dimensional PCs is 
warranted because of the inherent simplicity of the geom- 
etry and because a variety of extremely accurate mathe- 
matical tools are readily available which allow a detailed 
study of the problem without having to compromise ac- 
curacy. For instance, because the electric field and its 
first derivative are continuous across the interface, and 
because of the low dimensionality of the PC, the conver- 
gence problem that plagued band structure calculations 
for many 3D PCs 0>|g] is essentially non-existent for ID 
structures. Thus, we were able to use the old trusted 
plane wave (PW) method to find the band structure and 
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the DoS for supercell sizes not even imaginable in three- 
or even two-dimensional supercell calculations @>I3- O nc 
can obtain better than 0.1% convergence with as few as 
^30 plane waves per unit cell in the supercell. The trans- 
mission coefficient, too, can be calculated for nearly arbi- 
trary supercell sizes. Finally, one can calculate the band 
structure and the imaginary part of the wave vector using 
a semi-analytical approach for very large supercells. 

The ID PC is, in many ways, the "infinite square- 
well" problem of photonic crystals. It contains the es- 
sential features of its bigger cousins in two and three 
dimensions without the mathematical complexities and 
the accompanying numerical uncertainties [a, Q that can 
sometimes overshadow the essentials. For example, with 
3D face centered cubic structures, it becomes practically 
impossible, due to convergence problems, to increase the 
supercell size beyond 2x2x2 (or at most 3x3x3) con- 
ventional cubic unit cells which contain only 32 primitive 
cells per supercell (or 108), since typically at least ~ 1000 
terms per primitive cell are necessary to ensure sufficient 
convergence for inverse opal structures. It's not obvious 
from the start whether a randomness analysis with such 
small supercell sizes would yield results that are physi- 
cally meaningful. Artifacts due to the small supercell size 
are bound to be inextricably intertwined with the phys- 
ically significant bulk features of the imperfect PC. It's 
important to realize that, with the supercell method one 
still calculates the bands of an infinite perfect PC. The 
randomness is only within the supercell, but on a global 
scale, it's still a perfectly periodic structure! In order to 
be able to resolve the supercell artifacts from the physical 
features brought about by randomness, one needs to en- 
sure that the interaction between neighboring supercells, 
which, to a good approximation, is proportional to the 
surface area of the supercell, be small compared to the 
bulk properties of the imperfect PC, which can be taken 
to be proportional to the volume of the supercell. Hence, 
on purely dimensional grounds, one can argue that the 
surface to volume ratio of the supercell 1/L, where L 
is the linear size of the supercell, should be small com- 
pared to the typical length scale of the problem, namely 
the wavelength of the bandgap. We allowed the super- 
cell size N to vary from N — 2 to N « 9000, and one 
can clearly see the supercell artifacts gradually diminish- 
ing while the bulk features become more prominent in 
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the limit as N — > oo. On the other hand, ID structures 
can have features, such as a bandgap for any geometry 
and any refractive index contrast, that are certainly not 
shared by 2D or 3D PCs. 

The precise distribution of randomness in the geome- 
try of a PC would surely depend on the details of the 
specific manufacturing process. In the interest of sim- 
plicity, we chose the simplest distribution, the uniform 
distribution, in our study. As the unit cell, we chose a 
unit "supercell" that consisted of up to ~ 16000 unit 
cells. The thicknesses of the layers were perturbed by a 
given percent, by adding random numbers chosen from a 
uniform distribution. As the unperturbed structure, we 
chose the quarter- wave stack that has, for a given dielec- 
tric contrast, the largest relative gap between the first 
and the second bands, as can be seen in Fig. ^ I n what 
follows, we will consider this structure with a dielectric 
contrast of 13 as our perfect PC. For ID PCs, one fur- 




FIG. 1: The relative gap width vs the filling ratio / for a ID 
PC made of slabs of alternating dielectric constant of ei = 1 
and £2 = 13. The lowest gap has a maximum for / = 1 — 
\/ei/e2 = 0.72, which is the quarter- wave stack value. For 
this value of /, the even numbered gaps, the 2nd, 4th, etc, 
which are n general nonzero for an arbitrary value of / are all 
closed. 



ther has the luxury of calculating the bandgaps using an 
analytical method^. This approach also permits the cal- 
culation of the imaginary part of the wave vector in the 
forbidden gap region and allows a reliable assessment of 
the accuracy of the plane wave method for the problem 
at hand. 

We also investigated the transmission coefficient for a 
250 unit cell quarter-wave stack structure. The transmis- 
sion coefficient was calculated by simply matching the 
boundary conditions for the electric and the magnetic 
fields at each interface between the slabs in the multi- 
layer structure. 



II. DENSITY OF STATES CALCULATION 
WITH THE PW METHOD 

Maxwell's equations for waves propagating in the x 
direction in a medium with a dielectric constant e(x) that 
depends only on x, can be reduced to 



d 2 E 
dx 2 



1 . . d 2 E 



(1) 



where E is parallel to the slabs. With e(x) periodic along 
x with lattice constant a, and translationally invariant 
along y and z, 



e(x) 



a 



e(g)e i9X , with e(g) 



a Jo 



e(x)e 



: dx 



where g = m2n/a, is a reciprocal lattice vector with m 
0, =Fl, =f2, . . ., and E{x) can be written as 



E(x)=e ikx y, e ^v 



(2) 



where — 7r/a < k < it /a. For a given k, this yields an 
oo-dimensional generalized eigenproblem 



Q 2 E = —eE 

c z 



(3) 



or by multiplying both sides from the left by Qe , one 
obtains the ordinary eigenproblem 



(Qe^Q^QE) 



■{QE) 



(4) 



where Q = (k + g)S gg ', e gg > = e(g — g'), and e _1 is the 
inverse of the matrix e. For a given value of k, a trunca- 
tion of this oo-dimensional ordinary eigenvalue problem 
yields, by retaining only the g-vectors with \g\ < g ma x, 
the band structure u>j(k) and the modes Ejk(g). We 
choose a structure where the dielectric constant alter- 
nates between two values e\ and e 2 each with thickness 
dx and tfe, respectively. 

The choice of the lattice constant a = di + <i 2 is not 
unique. Although the choice a = d\ + di is the most 
obvious and the most convenient, the lattice constant can 
be chosen as any integer multiple of d\ + cfe, A = Na. 
With a choice for A with N > 1, and following the same 
formalism one can write 



e{x) 



and 



G 



e(G)e 



1 r A 

with, e(G) = - / e{x)e- tGx dx 
A Jo 



E{x) 



iKx E{x)e lGx 



where G — m(2ir/A), with m = 0,^1,^2,... and 
— it I A < K < it I A. Clearly, to get results with the same 
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FIG. 2: The band structure of a perfect ID PC with different choice s of supercell size N. The parameters of the structure are 
those of a quarter-wave stack, ei = 13, £2 = 1, and di/cfe = \J ei/ei- The points T and M of the "Brillouine zone" correspond 
to K = and K — tt/A respectively. When randomness is introduced, small gaps appear between each and every fold. 



level of accuracy as before, i.e. with N = 1, one would 
now need to include N times as many plane waves in 
the expansion, which simply increases the computational 
burden, both in terms of storage and computing time. 
The band structure for N = 1, 2, 3 and 250 are displayed 
in Fig for a perfect PC. i The folding of the bands in 
the first Brilloin zone for each N, makes the appearance 
of the bands rather different for each case, although the 
DoS and the eigenfunctions E would be independent of 
the choice of the supercell size. The frequency is plotted 
in units of 27r/a for all cases, so the frequency scale is 
not affected with the result that the bandgaps are at the 
same frequency, as would be expected. To calculate the 
DoS, we choose a uniform mesh in fc-space to calculate 
the bands and then choose a small frequency window, 
Aw, and count the number of modes whose frequencies 
fall within that window. 

We add random perturbations to the thicknesses of the 
layers in the supercell such that 



di, 2 = d? ; 



I _2y»(u-i 



(5) 



where d® 2 are the unperturbed values of the thicknesses 
of the layers, i.e. the quarter- wave stack values, u is 
a uniformly distributed random number in the interval 
(0, 1). We control the amount of disorder by varying the 
percent randomness parameter p between and 1 . p = 
corresponds to perfectly periodic structures, and p = 1 
corresponds to 100% fluctuation where d±, di can range 
between and twice their unperturbed values. When 
disorder is introduced, gaps appear between every fold 
for TV > 1. In Fig. EJwe plot the upper and lower limits 
for the lowest three bandgaps as a function of p, the 
percent randomness with a supercell size of N = 1024. 
Note that since for quarter- wave stack structures the even 
numbered gaps are closed, the bandgaps in this figure are 
in fact the first, third and the fifth bandgaps of a ID PC 
with arbitrary values for the layer thicknesses. The third 
gap centered at uia/2itc = 1.59 closes around p^ = 0.1, 
the second gap centered at wa/27rc = 0.96 closes around 




FIG. 3: The upper and lower band edges for the lowest three 
gaps calculated with a supercell of size N = 1024, as a func- 
tion of the disorder parameter p. 



P2 = 0.18, and the lowest gap centered at ioa/2-Kc = 0.32 
closes around p\ = 0.55. The ratios of the critical values 
of randomness p± : P2 ■ Pz agree well with the ratios 
of the corresponding center gap frequencies, C03 : LU2 : 
u>i. This can be understood using the simple argument 
that when the random fluctuations in the thicknesses of 
the layers become comparable to the wavelength of the 
gap center, the bandgap disappears since the destructive 
interference responsible for the existence of the forbidden 
band depends on the long range periodicity at that scale. 



A. Analytical Method 

As discussed in detail in RefQ, for n dielectric lay- 
ers with thicknesses d\, . . . , d n , with dieletric constants 
ei, . . . , e„, for a given w, one can obtain the transfer ma- 
trix, defined by 



Eq 
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(6) 
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FIG. 4: The imaginary part of the wave vector Ki and the band structure for selected values of the supercell size N for a 
randomness level of 10%. The band structure was calculated using the PW method and Ki was calculated with the analytical 
method. Note the slight shift of the bands due to PW convergence. As N grows, the second bandgap that lies between 
0.83 < uia/2nc < 1.07 in the perfect PC is more and more populated with transmission resonances, thereby narrowing the gap. 
For very large values of N, the bandgap appears to settle down to 0.89 < cua/2nc < 1.01. For iV = 4096, only Ki is shown, as 
the PW method isn't practicali for such a large supercell. The cusps in the Ki graph, which correspond to nearly flat bands, 
are still identifiable. 



where 



(7) 



with 
Di = 



1 1 



and Pi = 



e i^7Jujdi/c Q 

e ~i^eiu)di/c 
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Imposing the Bloch condition on the -E-field, one obtains, 



(8) 



Comparing with Eq. [fjl the eigenvalues of the transfer 
matrix are seen to be e lKA . Then, for t = \(Mu + 
M 22 )/2| < 1, K is real and is given by K = (l/A) cos -1 1, 
while for t > 1, if is complex with K r A = 7r and 
KiA = - ln(t - Vi 2 - 1) ■ 1 /Ki is the decay length of the 
evanescent mode, and is a measure of the strength of the 
bandgap. For finite PCs, it's desirable to have KiA 3> 1 
to have a significant drop in transmittance. The advan- 
tage of the exact method is that the supercell size N can 
be increased to values that are practically impossible us- 
ing the PW method. While with the PW method, using 
30 plane waves per unit cell of the supercell, the memory 
requirements scale as ~ (30iV) 2 , and the time require- 
ments scale as ~ (30iV) 3 , the exact method requires a 
very small amount of memory. The only disadvantage of 
the analytical method over the PW method is that, while 
in the PW method one chooses a real K and calculates 
the frequencies corresponding to that value of K, in the 
analytical method, one chooses the frequency w and cal- 
culates the real and imaginary parts of K, K r and Ki, 
corresponding to that value of lo. If the bands are nearly 
flat, as is the case for very large supercell sizes, then 



B1 



B 2 



E4 



2n 
V„ 

B 2 n 



B 2n+1 

En 



FIG. 5: Fields used in transmittance and analytic badgap 
calculations. 



one needs to sample the frequency interval of interest in 
very tiny increments in lo in order to "catch" a propa- 
gating mode. Thus the computation time can become 
very large. Also for large values of N, the transfer ma- 
trix M can have very large elements so one requires very 
high precision in order to calculate the transmission reso- 
nance frequencies. We used quadruple precision (128-bit) 
floating point variables and functions in the Intel Fortran 
compiler in order to be able to resolve the transmission 
resonances for supercell sizes up to N = 8192. For large 
values of N, even 128-bit precision is not sufficient, with 
the result that Ki cannot be made to completely vanish 
due to insufficient precision. Nevertheless, the propagat- 
ing modes appear as sharp cusps in the Ki vs u> graph 
which can easily be identified (Fig. 0J. For > 8000, 
one needs more than 128-bit precision to even see the 
cusps in (Fig. 0J. For larger values of N, we used Math- 
ematics for its arbitrary precision capabilities. However, 
compared to compiled code, Mathematica is slower by 
several orders of magnitude, so we had to stop at around 
N = 32000. To understand the bulk features of imperfect 
PCs using the supercell method, one would need a large 
supercell, in fact the larger the better. As the supercell 
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N 

FIG. 6: The dependence of lnT on Af for different values of 
the randomness parameter p for the center gap frequency of 
the first gap, wa/i2nc = 0.32. 



size is increased, small bandgaps begin to appear over re- 
gions that used to have propagating modes. One would 
normally expect the gaps of the perfect crystal to grad- 
ually shrink in size, rather than have more gaps, so this 
result seems somewhat puzzling at first sight. However, 
as the supercell size is increased, the statistical fluctua- 
tions decrease, and the pass bands become increasingly 
more densely populated. In Fig^we display the behavior 
of Ki as N is increased. As N becomes larger what used 
to be a photonic bandgap becomes more and more pop- 
ulated with transmission resonances, and the forbidden 
gap edges gradually approach each other, narrowing the 
gap. It's possible that as N — > oo, the whole bandgap 
region will be populated, albeit extremely sparsely, and 
instead of the bandgap, we will have a region where the 
DoS is extremely small — but non-zero nevertheless. We 
were able to increase up to N = 32768 and the bandgap 
was reduced as N became larger, although the decrease 
for very large N values was very small. To actually see 
the gap narrow even more would require an impractically 
large N. 

The number of transmission resonances in any fixed 
frequency interval Aw is proportional to N, so in the bulk 
limit with N — > oo, any wave packet with a small, but 
nonzero frequency spread Auj would contain many trans- 
mission resonances and thus would be partly transmitted 
and partly reflected. Clearly in the bandgap regions of 
the perfect PC, the density of these transmission reso- 
nances is extremely small, and these regions still appear 
to be bandgaps with a large, but still finite, supercell. 
Hence, it seems plausible to conclude that one cannot 
speak of a "true bandgap" for imperfect PCs, but only 
of large depressions in the DoS, which, in practice, would 
serve the same purpose as bona fide bandgaps. For in- 
stance, for a cavity made of an "impurity" embedded in 
a PC, localized cavity modes would eventually leak out 
through the PC "walls" of finite thickness, regardless of 
how perfect the PC walls are, bacause of the finite thick- 
ness of the walls. For such an application, what is im- 
portant is that the lifetime of the cavity mode be much 
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FIG. 7: (Upper left) lnT and Ki vs ui for a supercell of size 
N = 32 for a randomness level of 10%. (Upper right) A 
closeup for 0.8 < uja/2nc < 0.9 which contains the lower edge 
of the second gap 0.83 < ua/2-irc < 1.08 of the perfect PC. 
(Lower) Scatterplot of lnT vs Ki for the same structure. 



larger than the relevant time scale. Since the lifetime of 
the cavity mode is a function of the transmittance, for 
a given value of transmittance, one would simply need 
to use thicker walls as the random perturbations are in- 
creased. 



B. Transmittance 

Practical applications must necessarily use finite sized 
PCs, and for such structures, a quantity of more relevance 
is the transmittance. We calculate the transmittance by 
considering a PC of N unit cells, and each layer is per- 
turbed as described earlier. The transmission coefficient 
is calculated by imposing the boundary conditions for E 
D, B, and H at each layer boundary (Figure |SJ). This 
yields a set of 2n + 2 linear equations for the unknowns 
Ei, ... , E2n+2 ■ Setting the incident field Eq = 1, and 
assuming vacuum dielectric values for the incident and 
transmitted fields, eo = e n +i = lj onc obtains, R = Ef 
and T = E\ n + 2 for the reflection and transmission coeffi- 
cients, respectively. Alternatively, one could also obtain 
R and T from the transfer matrix as detailed in Ref. 
but with our approach, we can also obtain the Ei within 
each layer. For a given frequency in the gap region, the 
dependence of In T on the number of layers N is approxi- 
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FIG. 8: The DoS and InT vs frequency for various levels of randomness. InT (solid curve) is the average for an ensemble of 100 
random structures for each level of randomness. Also shown as dashed curves are InT ± a, where a is the standard deviation 
of In T for the ensemble used. 



mately linear for all values of the randomness parameter 
p, as is the case Q for the perfectly periodic finite crystal 
(Fig. EJ. However, depending on the level of randomness 
InT can change by many orders of magnitude. In prac- 
tice, this would mean that in order to obtain a given 
value of transmittance using an imperfect PC, one would 
now have to use a thicker PC. 

Although there is a strong relationship between the 
transmittance of a finite imperfect PC, and the modes of 
the superlattice formed by choosing the same exact finite 
PC as the unit supercell, this relationship is not perfect in 
the sense that the existence of a propagating mode does 
not necessarily imply a large value for T. Conversely, 
the existence of a bandgap does not necessarily imply a 
low value for the transmittance T. A close examination 
of the InT and Ki vs oj plots in Fig [7| will reveal that, 
on a large scale, the transmission has a dip where Ki is 
large, and when Ki is small the transmittance is nearly 
unity. However, a closer look at a finer scale in Fig. [7| 
one sees that InT can still be not as large as what one 
might expect from Ki. 

Looking at the plots of In T vs Ki Fig. one could ar- 
gue that Ki sets an upper limit for InT. This upper limit 
seems to be a linearly decreasing function of Ki. That 
the actual value of In T can fluctuate below a particular 
value for a given Ki can be understood if one consid- 
ers the coupling of the incident plane wave to the modes 
inside the PC. If, owing to the randomness, the propagat- 
ing mode happens to have a small amplitude near x n +i, 



it will not be transmitted effectively. Similar coupling 
problems were reported by Robertson el al [12| for 2D 
PCs. 



III. CONCLUSION 

We studied the behavior of the photonic band gap and 
the transmittance for an imperfect PC using the super- 
cell method combined with both the plane wave method 
and the analytical method. We also studied the trans- 
mittance for an imperfect finite ID photonic crystal and 
have shown that the results we obtain in all cases are 
consistent. The bandgaps of the perfect PC are replaced 
by a DoS that is extremely small to be detected even with 
extremely large supercells that consist of over 32000 unit 
cells. The higher frequency bandgaps disappear first with 
the lowest gap closing at around a randomness level of 
55%. 
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